#!/usr/bin/env python

import sys

import numpy
import histogram as H
import tables

from itop_reco import jf_map

hc = H.HistogramCollection()

hc.add(H.Histogram2D('e99r049jfocc', aranges=((0.5, 16.5), (-8.5, -0.5)), 
                     nbins=(16, 8), 
                     xlabel='jf-map x',
                     ylabel='jf-map y',
                     title='run 49 $q$ (all triggers)'
                     )
       )
hc.add(H.Histogram2D('e99r049jfocc1996', aranges=((0.5, 16.5), (-8.5, -0.5)), 
                     nbins=(16, 8), 
                     xlabel='jf-map x', 
                     ylabel='jf-map y',
                     title='run 49 $q$ (trigger 1996)')
       )

hc.add(H.Histogram2D('e99r101jfocc', aranges=((0.5, 16.5), (-8.5, -0.5)), 
                     nbins=(16, 8), 
                     xlabel='jf-map x',
                     ylabel='jf-map y',
                     title='run 101 $q$ (all triggers)'
                     )
       )
hc.add(H.Histogram2D('e99r101jfocc1996', aranges=((0.5, 16.5), (-8.5, -0.5)), 
                     nbins=(16, 8), 
                     xlabel='jf-map x', 
                     ylabel='jf-map y',
                     title='run 101 $q$ (trigger 1996)')
       )


bufsize = 1000
with tables.openFile("e99run0049m0.h5") as pt:
    integrated_q = numpy.zeros((4,4,8))
    for i in xrange(0, len(pt.root.adc), bufsize):        
        integrated_q += \
            pt.root.adc[i:i+bufsize].sum(axis=0).sum(axis=4).sum(axis=3)

    for col in xrange(4):
        for row in xrange(4):
            for ch in xrange(8):
                samples = pt.root.adc[i][col][row][ch]
                hc[0].fill(jf_map(row, col, ch), 
                           weight=integrated_q[col][row][ch])
    
                q1996 = pt.root.adc[1995][col][row][ch].sum()
                hc[1].fill(jf_map(row, col, ch), weight=q1996)


with tables.openFile("e99run0101m0.h5") as pt:
    integrated_q = numpy.zeros((4,4,8))
    for i in xrange(0, len(pt.root.adc), bufsize):        
        integrated_q += \
            pt.root.adc[i:i+bufsize].sum(axis=0).sum(axis=4).sum(axis=3)

    for col in xrange(4):
        for row in xrange(4):
            for ch in xrange(8):
                samples = pt.root.adc[i][col][row][ch]
                hc[2].fill(jf_map(row, col, ch), 
                           weight=integrated_q[col][row][ch])

                q1996 = pt.root.adc[1995][col][row][ch].sum()
                hc[3].fill(jf_map(row, col, ch), weight=q1996)

                
print

H.dump(hc, "jfocc.hpt")


            
